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A HYBRID MATHEMATICAL MODEL OF COLLECTIVE MOTION 
UNDER ALIGNMENT AND CHEMOTAXIS 

E210 D1 C’OSTANZO* AND ROBERTO .NATALINI^ 


Austkact. In thi& pc^jct wc propose an.i £tiidy ft hybrid discrett in eontit&uout nuillicniftlirnl 
modc'l of coUcctive uiouon under alig;iun<'Qt and cbcinotaxift udeci. .Siartlug frorn liiu paper 
by |Di CosiangQ ct al |J20^, in whkii the Cucker-Smalo model |JCiJc|«|ranrlSro^ 200^ was 
coupled witli ocLor coll muchmiiMiift. lo describe the cell iJiig;rftilun and sulT-orgauisftiiun in 
iBe zebradah lateral line primordiimi, we inlrcKluec a simplUied model in whidi the coupling 
between an aligmnent and cliemotuxu) uiccbanksm acts on a system of interacting particles. In 
particular we rely on a liybiid tlcscription in which the agents arc discrete entities, while the 
dieinoaitractant is eunaiderod as a continuous signal. The ptupused model is then studied both 
from an analytical and a numerical point of view. From Uic analytic ]>oint of view we pro\‘o, 
globally in time, existence and uniqueness of the solution. Then, the asymptotic behaviour of the 
system is investigated. Through a suitable Lyapunov functional we sliow' tiiat, tor thu linearised 
version of our system, tor t -* ^oo, the migiating aggregate exponentially converge to a state 
in which all the particles have a same position with zero velocity. Finally, from thu nuinuriual 
point of view, sfimo meaningful dynamical tc«ts are proposed to simulate the behaviour of the 
system, alwi m comparison witii Uio amJvtical dudingS' 


1. Introduction 


A collective tnotiou is a form of collective behaviour in which the individual unit’s action is 
strongly dominated by tlio influence of other units, so that its motion results very different from 
how it would be if it was alone. Collective motion has been extensively studied hi lecent year 
in a great vaiiety of systems, from nou-hviug systems, such as nematic fluids, nano swimmei's, 
simple robots, to living systems, such as bacteria colonies, cell aggregates, swarms of insects, 
flocks of birds, schools of fish, swarms of insects, herds of mammals, crowds of humans. Many 
references can be found, for example, in the review by [Vlcsek and Zafeiris ( ]2012[ |. 

Ditfeieiit reasons of aggregations in a given system have been proposed Coordinated 
motion of cells results in making a biological process more efficient (e.g. in embryogenesis, 
wound healing, iunnime response, etc.), while in the case of tumoiu* cell invasion it appeals to 
speed up the progression of the decease jMehes and Vicselc|[2Q14[ i. For the animals entities it is 
ob6er%'od that a group can more efficiently explore surrounding environments, to enhance foraging 
capability and detection of piedatoi's ^Pitcher et al^[I08^; being in a group oilers piotection 
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against attacks by predators jloanuou et al 2008[| aud mcreases the lucomotioii cfficipucy j]Fish 


l9j^. 


The* research iu the field of collective motion modelling is of gioat interest for the 
applications in many tields of the real life, from biotuedical fields e.g, iu relation to the collective 


20061 IBebnunte et all 120081 lArboleda-Estiidillo et all 12010 


bepUlvoda et all 2 Uiat |.loie et all 

57115' 

Loliu et all |2ui3t |Ui Costauzo ct all | 2 (Ubl and the 

review's by 

Hatziklrou aud Deutsch 

2007 

and 

Mehes and Viceek 

2014J, to eocio^ouomic field 


2010 Pareschi and Tuscaui 2Q14[| to the problems of pedestrian 


flows OHelbing et alH1997| |Plccoli and Tosin[|20Q9| [Faria et al[|20f0| |Criatiani et all |2Q1Q[ |Bruno 


etall 12011 1 IMoussaj’d et alj 1 20111 ICristiani et alj 1 20141 |2QI5>. 


Basic mathematical models about collectIw motion are substantially based on one or 
more of the following steering behaviours of the units: alignments iu or<ler to mot'e in the same 
diiection of average heading direction of neighbuuiing units ( JVicsek et al[ 1095 Cucker and 
Smale 2007 and references therein); separutiuft, In order to avoid crowdiug!''a3I^ cahtaioti. to 
remain close to the average position of neighlx^ms ( jP^Orsogna et a4 2006[ Strbmbom[ 2Qll[ f. 
Sometimes three-zone models combine together alignment, attraction and repulsion effects on 
three uou-owrlappiug regioiLSs w*hich take into account the sensory capabilities of the individual 
JAokii 119821 |Cou 2 u> et~al |2002t [Gregoire et alj |2003i Gregoirc and Chatel |2Q04t l. For example, 
the vibu^ hold of a bird does not extend behind its body^ Ashes can acc^ompany visual signals 
with those coming from them lateral line, while cells can feel stlmuh around them. This leads to 
introduce the concept of the cone of vision j ]Huth and WLssel 1992 Couziu et al[|2Q02] Hemelrijk 
aud Hildenbrajidtl 2Q08||. In general the cone of interaction is not only aulmal-dependont, but it 


can also vary depending on the type of motion, environmental conditions, presence of predators, 
and aim of the displacement ^Cristlanl et al 2Q1I[ |. Moreover, when each unit interacts only 
with those particles which aie closer than a predefined distance, that is a lixed range of influence, 
we speak about “metric'* interaction. On the other hand, I'ccent studies of starling flocks, on 
the other hand, have sliown that each bird modifles its position relative to the mates directly 
surrounding it (typically six or seven Imllviduals), no matter how close or how far away those 
fluimal.n In this type of behaviour, referred as “topological' interaction I jBallerini et al 2Q08[ |, 
the radius of perception is adjusted by each individual, in such v^^ay that the neighbourhood of 


interaction encom 

lasses a predefined mass of other individuals felt comfortable to interact with 

It |)C.'nstiaui et al 

2010 

lUristiaui et all 121111 and references therein). 

VVltn rcgarc 

to t 

»e alignment models a widespread model is represented by the Cucker- 


The original model was proposed to describe the 
dynamics in a flock of birds, but its appUcability is to general phenomena where autonomous 
ageuts reach a consensus, e.g, animal hertling, emergence of common languages in primitive 
ijCouzin et al 


societies, etc, 
seminal paper by Vlcsek et al 


aud Jadbabaie et al 


2005 ( |. The starting point of thi^ model Is represented by the 
j]1995^ au<l prevdouH analytic studies cau be fbimd aLso in 


Anyway, the Clucker-Smalc paper is widely known to have 
established an aualytlcal Txact result on the convergence to the same velocity In a group of 
interacting ageuts through an alignment effect. 

The hypothesis of the Ciickei-Smale iiiodol is that the force acting on every particle (bud) 
Is a weightc<l average of the ditferenccs of its velocity with those of the other particles (binls). In 
particular, for .V € N particles, the proposed model, in c^outinuous time t € R>o« can be written 
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ill the form 


1 ^ 




-r(Vj-V, 


,^i(a2 + \\X,-Xjfy' ' " (»•!) 

. X, = V.. 

where X,s V, € R'‘, for i = I. N. are the position mid the velocity of the i-th particle, fti, 

ci 2 . c are positive constants and ||-|| is the Euclidean norm in R". In particular c captures 
the rate of decay of the influence between agents in the hock as they separate in space, and it 
is a fundamental paianieter in the time-asymptotic liehaviour of system (see below). In 
particular, under suitable initial data and parameters the model implies that the state of the 
group converges to one in which all particles nio%'e with the same velocity, said ^odjWn^ sfafc, 
The main convergence result proved in Cucker and Smale | )2QQ7[ f has been unproved by |Ha| 
using an explicit Lyapunov functional approaHi. First, as in this last reference, 


and Liu 


we can give the dofmitiou of tiine'OJiymplolic flocking as follows; 

Definition 1. Given syittem (jOJ. Ut X('a/ = ^ X< and VcM = poslfijn 

and velocity of Vie cctifrc of mass. System j ]l.l| i has a time’asymptotic flocking if and only if 
(X,, V,). i = 1.,,.. satisfy the two conditions: 

i) the velociiy fluctuations go to zero time’CsymptoliaiUy (velocity alignment): 

,v 


_Un,^^||V.{O-Vp„(t)||= = 0; 

|sl 

3} the position fluctuations arc uniformly bounded m time t (forming a group): 

yjX.{t]-Xr„{t)\\‘<+cc. 


sup 

0<t<+oo 


1=1 


( 1 . 2 ) 


1.3) 


Notice that the squoie root of the quantities under the limit and supremuin operations 
in ( ]1.2^ i |1.3[ ) is proportional to the standaid deviations of Vj(t) and X<(£) aiound the centre 
of mass system. Cucker and Smale) j ]2D07t and [Ha and l1u| ( ]2Q00[ ) proved that for s £ |0,1 j2\ 
occurs a global unconditional flocking of system ^l.l| |, as stated in Dehnitioii regardless of 
initial (x^nligurationH. while for a £ (1/2,+cc.) there is conditional flocking, that is only some 
parameters and initial data lead to a hocking state, but hi general the dispersion of the hock 
may occiu'. 

Cucker and Smale-like models have been widely employwl iu collective dyuauiics and 


several applications can be found, from the biological held 

Szabo et aJ i 

lelmonte ct al 

2(11)8 

Arbolcda-Lstuclillo et all |2U1UI 

hicpUlvoda et all 


Ui LJostauzo et al 

to the 


extensions have been proposed. For example, [Cucker and Dong )]2010[ 2011|l added a repelling 


force betw<*en particles in equation proving au analogous convergence theorem eusuriug, 


on certahi conditions, flocking behaviour, iu addition without collislou between particles. IHa and 


Levy (2Q09[) apphed system jjl.ip to describe the motion of pbototactic bacteria (i,e., bacteria 


that mo\'e towards light), adding a force that mo\'c excited bacteria towards the light source iu 
equation (jOJi, aud an addltloual equation for the excitatiou level of each particle. Also in this 
case, uudei particular conditions, it Is proved that the asymptotic velocity of the particles tends 
to an identical terminal velocity. Other models have hitroduced noise jCucker and Mordecki 
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201)8[|. Stochastic equ ations j ]Hd cl al| 2000( 1. leader mdlvldual^ with a preferred beadiiig dlrectiou 
JCucker fuid Huepe[ 12006( 1. kinetic eqxiatioiiti dCarrillo ct al |201Q(|, 

In this paper, starting from the model proposed by |Di Costanzo et al ^2015( [ In relation 
to tlie morphogenesis In the zobrahsh lateral line priuioidium. In which the authors coupled the 
Cuckor-Smale model with other cell mechanisms (eheiuotaxis, attraction-repulsion, dainping 
effects), T.'e study, both from an analytical and a numerical point of view, a simplified inodel in 
which the coupling between an alignment and chemotaxis effect acts in a system of interacting 
particles. Our description is hybrid; discrete foi the hiteracthig agents, and continuous for 
the chemotactic signal, FVom the analytic point of view, in the two-dimensional case witli N 
particles, we prove, globally in time, existence snd uniqueness of solution of the model. Then, 
the time-asymptotic behaviour of the inodel is m%‘cstigated on a Uiiearized version of the system. 
In particular, through a suitable Lj’apuno'.' functional, is shown that position and velocity of 
ail the particles go exponentially to those of their centre of mass. Moreover, the velocity of 
the centre of mass tends time-asympCotieally to zero. From a munerieal pohit of view, on a 
bounded domain with periodic boundaiy eouditious, some 2D dynamical tests ate proposed to 
simulate the behaviour of the full uouhiicar model, au<l to compare the uumerical results with 
those aualytical, 

From our findiugs, we observe that while in the Ciuker-Smale model two conditions may 
occur', conditional and unconditional hocking. In our model the hocking behaviour, as given In 
Definition!^ is always ensured. Moreover, we haw the stronger conditions that all the particles 
C'onwrge In a single position, i.e. their centre of mas.s, and that the velocity of the centre of 
mass tends to zero. 

Among the various numerical simulations, we discuss the competition between ali gn ment 
aud chemotactic effects varying the parameters of the system. We find a decx'case in the rate of 
convergence of the particles when the sti'cugth of the aligtuucut term increases with respect to 
the diemotaxis. On the other haud, an increase in the rate of convergence can be found with 
the same parameters aud with a greater number of iuteractlug particles. Inspired then from 
the aforesaid paper ou the zebrahsh lateral liuc, we cousider also the case of two kinds of cyclist 
leader cells ciinespondiug to the sources of the diemotactlc signal, and follower cells that do uot 
produce any chemical signal, both subject to the alignment aud chemotaxis force. lu this ease 
wc observe the convergence of the group toward the soui'ces of the chemoattractant. Finally, we 
c'ousidcr the collective motion under chemotaxis. ueglecting the alignment effect, Our uumerical 
tests, in this regard, show the alkscnce of convergeuce, and an oscillatiiig motion of the particles 
around tlieir centre of mass. This suggests that, in such model, the only chemotactic effect is 
unable to reproduce biological phenomena involving stationary aggregates. 

The paper is organized as follows. In Sertioii we design our hybrid mathematical 
model of collective motion. Section deals with the analytical study of the local existence and 
uniqueness of the solution. Then, In Section 2] this result is extended and the global existence 
aud uniqueireas is pro\'ed. Section is devoted to im'cstigate the asymptotic behaviour of the 
model around the equilibrium couffguiations. In Section l^wc show and discuss some meaningful 
dynanucal simulations of the model. Finally. .Section [^includes the conclusions and possible 
future pei'spectives. 
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2. The basic; mathematical model 


Staitiiig from model we suppose that the force acting on each particle is given 

by an alignment term, proportional to the dlifereuces of velocity with the other particles and 
weighed on the distances« and by a chemotactic attraction towards higher (^ncentratlon of a 
chemical signal /(x, £). supposed produced by the paiticles themselves. Typically, this last foi'ce 
is proportional to the gradient of the concentration V/ {see Eisenbach and Lengeler| 2Q04| for 


biological baihgrounds, and Murrayl 2003 Perth;4inet 2007 for some mathematical references). 


In our hybrid description, while particles are consider^ discrete entities, the signal / Is supposed 
to be continuous and its rate of change in time la equal to a diffasion term, a source term, and a 
molecular degradation. If each particles Is endowed of a radius fl. describing its ciieular shape, 
the source of signal can be given by a characteristic function on a ball of radius R ceutiad on 
each particle. 

To summarize our hypotheses, we w'ritc the following system: 

s 


V 




1 + 


|X.-X,| 




X, = V,. 


( 2 . 1 ) 


S' 


dtf = DAf + C ^ Xb(X^ R) - V/- 
where ^,tT.'y.D.^,T}.R are positive corrstanCs, and 


_ r 1. if X € 

^ Q otherwise: 


with 


B(X,.i?) := |x:j|x-X,|| <i?}. 

Initial data are given by initial position and velocity for each particle: 

X(0) = x„. 

V(0) = Vo, 

with X = (Xi...,, X,v), V = (Vi,,.,, Vat), and by the Initial concentration of signal, that we 
assirme 


/(x,Q) :=/d = 0. (2.2) 

We note that equation j |2.i^ 3 can be analytically solved making the classical exponential 
trartsformation; 


S' 

dtu = DAu + ^ e '" ^ \B(X^.ft). 

•-J 


with ^(x, £) solution of 
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aud u{x.Q) = /l). Now, if 


, . I UiNf 

is the fuudamoutal solutlou of the heat equatiou on R^‘, using the Duhamol's principle we have 


/(X, 0 = (r(x, t) * fo)e-^ + r(x, t-r)* VB(X,(r),R) j 

^ r< / 

= 4 ^ e y^ r(x - X, i - r) (t-) j?) t^xdr 


dr 


(2.3) 


where » is the c^onvolution operation In the variable x. and /o = 0 for our Initial C'ouditiou j ]2.2[ f. 

In the following, for analytical and nutnerlc^ simplicity, we will discuss the case of N 
particles in R^, First, if x = (xi, 2 : 2 ) atJd x = (xi.xs), from j]2.3tl we can write 


^ rl 


f{xi 


^>0 y/BfX, 


g-rfft-T).- Ht-r\D - 


frl-ft) 




dxi dx'i dr. 


(2.4) 


With a simple change of variables It Is easy to see that the integrand In | |2.4p is with I'espect 
to r in (O.t). In fact, if we define 


I'i - xi . xa - X2 

. . . =, xa = 


V4(f-T)D 


V4Ct - t)D’ 


(2.5) 


and 


:= (-ii%/4(7^^r)D + i2,-xa\/4{t^^T)D + X2j » 

we have the Jacobian 


det [d^/'dy^ = 

In the new variables < \2A\ becomes 

A' ft 


-^/4(r^T)D 0 

0 -y/i{t-T)D 


= i{t-r)D. 


/(Xl,X2.f) =C 




4?r(C — t)D 


si//. 


-4(l — r}D<ixj (1x2 dr 


(2.G) 


where B Ls ubtabied from B under transfoiiuatlou In i )2.()[ l the integrand is L* in (0, f), 

with respect to r. 
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Now, if V/ = ^2/) ^ chemotactic gradieuts aud X, = (X,isX»2)> V» = 

are renpectively the position and the velocity of the i-th cell, we obtain 


A* 








,-*7(t-r)g -4(e-T>C 

4 jT(t — r}D 




.£ V / _ / e~ 

^Jo 4ir(t - r)i?/A:,2(r)-ft 


ri(r) -(i2-(r))* - xi) 

X / e we-'j _1--^<iridx2dr 

yY,.(r)-Vft^-(x2-X,*(T))2 4(t - tP 


^ ‘f _<£2;_£2if 

£ 4ce-')o 


47r(f - T)i) }x,,{r)-R 




X e 


tttttttt 


— e 




cti'^ dr 


iL r «-*?<*-’•) y+ft (^-x,a(r>-*2>a / 

fV' / _ / e-- \fi - 4ii~r)D - 

4ir(t-r)Oy_ft 


-e j dx 2 dr. 


where iu the last step we have used the change of variables 


X2 = X2- Xj2{r). 


{VtlsV,2) 


(2.7) 


Similarly we can proceed for dxjf{^i^X2 t)- 
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Fiually. substitutiug i )2.7[ f iuto ( j2.ip i wc caii summarize, for i = 1,_ N. the following 

system: 






t2 


a ^ 

= -T 

N ^ 


j-i (l +-TP 




x,-x, 


e 4(t.T>0 




fi) - ,\j , (r) - ^/^Tj3y 
-- 


- 1 di2dT, 


a ^ 

=-y 

AT ^ 


r< g-'iU-’-) ^ r+« 

■'V,. 4 ^(( - 7)D ^ 


(Vj2 - v:,) 


( 2 . 8 ) 


g 4(t-V)0 


(>\210 - 

- Ht~r)D - 


(aj <«). Xji (T).’ 



—e 

Xa 

• 

II 

X.2 

= V^2. 


- dr. dr 


3. Local existence and uniqueness of the solution 

In this section, using a fixed polut argumeutatiou, we prove for ( ]2.8^ the local existertce 
and miiqueuess of the solutiou, lu the next section, theu, this result will l)e extended to a global 
one. 

First, let y = (Vi.Vjs , Xi,_X,v) the solution vector, and let 


q:= (qi.qj), 
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with 


and th^u 


with 


Finally, if 


9i.»i -= 


a ^ 

-T 

AT ^ 


j=i \1 + 




P 




(^2-^2). 


92,»1 - 
92,^ -* 


v;i, 

^2- 


P-= (Pl P2)i 










—e 


(x., (0-A^,W + ^R®-i|^ 

^ ] di'i. 




• U-r)D 




—« 


(v^ t j) - ^jS If) ' 

- 1 


P2,»l =P2.i2 =0. 


(3-1) 


(3-2) 

(3.3) 

(3-4) 


(3.5) 


(3-0) 


~ 4jr(f-r)i5‘ 

system ( ]2.S^ can be written as 

y = q(y)+ f C{t - T)p{t -T.y{t),y{T))dT. 

Jo 

lutegratlug from 0 to t we have 

y = yn+ f q(y(r))«fr+ f f C(s - r)p(s - r,y(s).y(r))<fTds, 
with yi) = y(0). Then, interchanging the order of integration in the Kecond integral, we hax’e 


y = yD + jf'[q(y(r)) + |^'c(s 


-r)p(.9-T.y(s),y(r))ds 


dr, 


(3.7) 
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or 


y = yo 


f [q(y(r)) + h(«,T,y(T))| dr, 

Jo 


(3.8) 


with 


h(t.r-y(r)) ;= £ 


C{dJ-r)p(A-T,y(s).y(T})Js. 


For A discussion of such type of equations see^ for example, |Burton| (l2Q05t. IKhallll (]2Q02t. and 
also [Wazwa2| { |2011[ |. |Laksbmikantham and Ramil i ]1905p . 

Now, lot > 0, We consider the set 


S' = {{tT,s,y) :0<r<a<t<a, ||y(f) -yoH < t>}. 
Since q(y) is continuous on S, we can define 

A/, =max||q(y)||. 


(3.!)) 


Then, to prove that h(?,r,y{r)) is coutlmiuus in S, firstly we prove that C{s — T)p{s — 
T,y(s).y{r}) in L^{r,t) with respect to variable a. It is enough to deiuuustrate the iiitegrability 
aioiind A = T. Starting from we consider — y(A).y(T}). and the change of v'aiiables 

~ ~ (310) 

Since the quantity in brackets is bouiiiled by 1. we have 


|pi,ii(s-T,y(s) 


S 4i.-rfO rij-2 


.y(r))|<X^ / 




“TTTTTrr 


v 


j^lJ-oc 

= A7\/4(s —t}D, 


with 


^f := N j e ‘ dt. 
J-oa 

The same holds for pi ,2 in j]3,fip. so we write 

aud 


|P(» - r,y(»),y(T))|| < ^Iym^/^s-T)D. 


|C(a-r)| ||p{a-r,y(s),y{T))||< 


2ffV(*-’-)D 


(3.11) 

(3.12) 

(3.13) 
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li 


Now C{3 — r)p(fl — T, y(s).y(T)) is (^otmuous in y and, from ( 3.13[ |, it is L^{T.t) with respect 
to variable s. so h(?. r. y(T)) is cootuiuous in 5. and we can define 


A/2 = mp||h(t,r,y(T))||. 


(3.14) 


To prove local existence au<l uniqueness, we want to obtdu a Lipschltz condition In S for the 
functions q and p with respect to the variable y. First, because of q is C* on S. the Jacobiau 
matrix (c)q/dy] is boimded on S uniforuily in r, so q satisfies the Lipschltz condition 


l|q(yi) - q(y2)ll < llyi - ysll. 


with Li positive coustant and (Lr, s, yi), (1 r, s,y 2 ) 8 S. 
Thou, from )]3,rj[( we define 


(3.15) 


- r,y(fl), y(T), xi) t* e 


f «).Xj ,(T)- y/ffi 






— e 


4l«-T]n 


(3.16) 


For the Jacobiau matrix of hold the following equalities 


dVki dV^2 


= 0 




(J) - Xj , (r) ' 




-- 




- k = i', 


ijk 


—e 


Xul 


2[»-T)P 


X,i (J)- X ,, (r). , 

- irr-^ -^X«(3)-X,i(r)+V^735T 

‘ ^ - i{s-r)D - 


. * 




= 0 . 


dXkz 

Since these quantities are boimded in S uniformly in ii- ■'*> we can write 


|pi,tiji(fi-T,yi(s),yi(r),i 2 ) -pi,uji(s-r,y 2 (^).y 2 (r),i 2 )| 


<^2(l|yi(A) -y2(s)||+ ||yi(T) -y2(T}||), 


(3-17) 
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where ij is a pasitUti constaut, and (r. s. yi), (T,3,y2) € S. FVom this ^*6 deduce 
- T-yi{.s),yi(r)) - ny 2 (a), y 2 (r))| 


t/: 






— e 


TUlsTry 


X - T,yi(.s).yi(T),f2} -Pi.tiji(a -T.y2{a),y2(T).i2)| di2 

=g 3 Q \/4£>{s - T)l 2 (llyi(ii) -y 2 (ii}|| + ||yi(r) -y 2 (r)||) 


N / .*'•3 






(tHR 


V ./ITTTTTB 


,a 

V^C—r)D 


r+cc 

< %/4D(a-rV2(||yi(a) -y 2 (R)|| + ||yi{r)-y';(r)||)2^ / e"* d 2 

v'4D(s-T)2A-n2{||yi(.^) - y2(^)|| + ||yi(T) - y 2 {r)||), 

where we have used the change of vanal)les | |3.10^ . Moreo\'er, are the variables 

belo ngin g respectively to the vectors yi. y 2 i and (f, r, Ssyi). (Is t. a, y 2 ) € S. 

The same can be doue for pi,tit so p satisfy the following condition in y on S': 

I|p(fl-^>yi(s).yi{r))-p(s-r,y 2 (s).y 2 (r))|| < ^/4D(j - 7)1^ (||yi(<^) - y 2 (s)ll 

+ ||yi(T)-y2(T}||). (3.18) 

with L 2 a suitable positive constant that Incorporates previous constants, and (r, s,yi}, 

(r,s.y2) € 5- 

Now, we fix 

l> 1 


T = min n, 


Afi +Af2'Li + 2l2A/ 


(3-19) 


with il/i- A /2 given by ( |3^ ( ]3.Up > and 

A;:*^^*|C(2)|>/4^d2. 


(3.20) 


Then we prove the following 

Theorem 1. E^udliofi ^3.8[ i has a umtjue solutwn on (Qs7^. T w defiritui in ) ]3.19^ . 

Proof. Wo consider the functional space 

S = {yeC“([ 0 ,r]):!|y-yo||c 70 < 6 }, 

where 

||y-z||c«:= sup ||y(I)-z{I)||, 

0<t<T' 

aud tte define the functional A . ^ > B as 


A(y}(I):=ya+ f [q(y(r)) + h(f,r,y(T))] dr. 
Jo 
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To spc that A : 6 —^ S notice that y continuous implies A(y) coutluuous, becaasc q aud h are 
c'OQtinuous, and that 


l|A(y)-yollco = sup ||A(y)(0 - yi)|| 

D<t<T 

< sup f (||q(y{r))|| + ||h(t,r,y(r))||) dr 

D<t<r JO 

< {Ml + M2)T 

<b, 


where we have used aud. lu the last inequality, 13.19To see that A is a coutractiou 

mapping, uotlce that if yi and y 2 € B theu 


mappmg, uotice tnat ii yi auu y 2 
l|A(yi) - A(y 2 )||c^ = sup || A(yx)(f) - A(y2)(0|| 


0<t<T 


< sup 
o<t<rJo 


[ lllq(yi(T))-q(y2(T: 

Jo 

+ j |C(s-r)|||p(*- r,yi(.<),yi(r)) -p(.v-T,y 2 (s)-y 2 (r))||)dsdT 

<(™™ sup [ [Li\\yi{r)-y 2 {T)\\+L 2 f |C(.s - r)j ^4D(s - r) 
^- - 0<t<TJo Jt 

X (llyi(5) -y2(-9)ll + llyi(T)-y2(r)||)|d.4dr 

< LiT llyi - yallco + 2I2 |[yi - y^Wco 


sup / / |C(s - r)j Jw{s - T)d.<idT 
<t<7Vo Jt 

LiT + 2L2 sup f f jC( 2 )| \/4Dld2dT I ||yi-yall^ 
I o<i<tJo Jo } 


< (^LiT + 2L2 I*" \C{^)\y/^:d-.dTj\\y,-y2\\(y 

= {L} + 2L-iM)T ||yi — yallco • 

Fvom p.inp the coustant (Li + 2L2M)T e (0,1}. The Banax'h-Caccioppoli fixed-polut theorem 
completes the proof. 


4. Global existence and i;niqi;eness of the solution 

To obtdu global cxlstcuce and uniqueness for ( |3.7p we will use a principle of continuation 
of solutions. We will pro%'e that bounded solutions can be continued to t = + 00 . The following 
general result, adapted for equation i |,1.7^ . provides a condition for the continuation of solutions. 

Proposition 1. L£ly{t) bt a fiolution of ( ],3.7^ on a intei'vai [0, T), ifViere w a P usilh 

||y - yoll [0^7"), then theie i‘i at >T 6'iich that y{t) can be continued to\0,f]. 

Proof. We show that lim^_* 7 '- y(t) exists, so we can apphed Theoi'emj^ star ting at t = T, 
aud this completes the proof. 
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Lpt tn bp a moiiotouic iucrpaaiug sequencp with limit 7". »uid let 

U = {(t.T,.^.y):0<7<s<t<T. ||y-yo||<P}. 
We prove tliat {y(Ch)} la a Caucliy :>eqiieiice. If tm > C»t. from fl3.7li we have 


||y{(m) - y((»)|| = 1^ [q(y(r)) + ^ - r)p(s - T,y(s), y(T)] 


ds 


dr 


C{3 - T)p{s - T.y{3).y{T))dji 


< I lly^ C(.H-r)p{a-7,y{a),y{T))ds 

- C(s-r)p{a-r.y(a),y(T))d« dr 

+ 1^ [q(y(r)) + j C{s - T)p{3 - T,y{3).y{T))d3 

a-r)p{a-r.y(a),y(T))d«| dr 

+1!^ / 

+ 1/ q(y(’'))rf’' 

lUtn 

< / / l|C(5-T)p(a-r,y(s),y(T))|| darfr 

Ju 

+ / llC{a-r)p(s-r,y(s),y(T))|| dadr 

rlH> 

+1 llq(y( 


dr 


dr 


dr. 


lu the last inequality the tliiid integral touds to zero as n, in —^ +oo, becaus<' q is bounded on 
U and fmtffi -* Also the first two integraLs go to zero as n, m —^ +cc, because of ( 3«13[ f. The 
proof is completed. 


Now, from Proposition [T] we obtaiir the following 
Theorem 2. Eqaalion has a unique global solution for ailt >(f. 
Proof. Fu'Kt. equations ( ]3.1^ ( ]3.4^ imply 

<23||y||, |9U2| <2/?||y||, 

< \\y\l \Q2.a\ < llylL 


so that 


||q||< V2A'(4-?’ + l)||y||. 
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Then 1)3.12^ establishes 


p(s - r,y(a),y(T))|| < My/WyjAD{3-T), 


so from we have 


< llyol 


V2iV(4/P + l) rilyWII 


dr 


mV^J' J' |C(i - T)\^/W{7^dsdr 


||yo|| + V2iV(4/P + l) /'llyWIIrfr + A/v^ f f '|C(s)| 


y/ADididr 


||yo|| + V2^(4^’ + l) f \\y[r)\\dT + M^ t \C{z)\ViD 

Jo Jo Jo 

\\yQ\\ + ^/2N{A^TT) J ||y(r)||rfr + A/A7V2]Vt, 


4D5 dz dr 


where ^*6 ha%'e set r = a — r. aud A/ is giveu by \ S.' 
Now, for each 0 < i < T. 


I|y|| < (llyoll +A/A7Vw7r) + V 2 jV(ifl^^'||y(r)||dT, 


so that 


||y|| < (||yn||+A/A/V^r) 


by tbe GtoiiwaJl’s inequality. Since soliitiuu remains bounded, for Proposition [TJ it can be 
continued to all (0.+cc). 
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5. Asymptotic properties on the linearized system 

In this section ^*6 prove Komo asymptotic properties ou the liuearized form of system 
that for convenience we recall here; 


Vi: 


N 




- , nx.-x,i) 

1 + -JJT 


;=i{ 

rt r+i 

/ C{t-r)Y, / 

Jo fTi J-i* 




^ r+R (x.ain-x,2ir)-ia)2 / -a,. 






—a 


(.V„ U) - STj, Ir) 

) diqdr. 




= -E 


7^{Vj2-V,2) 


j=l i J H- jp - 


Jo ^y-R 


4tf-r)D 




—a 


(,Yc2 (0 - Xj2 [T) * ^ 


y^ii — Ki, 


(5.1) 



= v;2, 


with 


C{t - T) 


7Ce 




4jr(i — r)D 

We axe interested in the equilibrium points that satisfy the condition; 

X,(?) = Xeq(f), Vl, Vt; 


(5.2) 


v,(n = 0. 


V:; 


^ X, (f) = = constant. Vi, 


(5-3) 


Equation means that all partirlos aie in a same position for all times. Now, to make a 
first-order appimimation of )|2.8[t. we cuusider the following Thy lor expansions around points 
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F, (Xj - X„ - v,i) ■.=-^ -_ v.i) 

" j 

= Vi:-Vu + Pi{X,-X.Xji-X\,), 
F2{t-T.Xa{t)-Xj2{T).i2)-.=£ -^TTTtTB- 

4 

ITTTTTJ + - r, X,2{t) - Xj2{r), ig). 

(aT.iW 

F, {t - t,XM - Xj:{Tli2) :=«-- 


=e 


Mt r)D 




• («• ')D 


2{t-T)D 

+ P3(t- ’•.X.i(t) - -V3l(T),i2), 


(X.,(0--V„(T)) 


(5.4) 

(5.5) 
(5.0) 

(5.7) 

(5.8) 


(5,!)) 


where the functions aud ^ contain the uunJiueai terms, while ^2 coutains the lineai' and the 
uunlmear terms. Similarly we can treat equation 

From 1 )5.jp 1 )5.flp , we liueari 2 e equation ( j.Vip i in the form 

= 

iv 

x'^[Xn{t) - Xji{t) + X,i{t) - Xji{7)) <1x2^7 

j=> 

jV 

x^(A*,j(7)-.V,,(T))dr 


'H?-xldl2 


N 


/ C(f-r)y (A-nU)-AV(r))dr. 

j=i j=i 


with 




C{t - T) 


2(f-r)D 


(5.10) 
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Similarly it can bp done for Via. Finally, we obtain 

JV .» JV 


^‘ = I - K.(t)) - E iX'iW - >^Ar)) dr. 


N 


J=1 

N 


^*2 =4E(V;2(i)-K2(0)- / C(t-r)^(X„(t)-X,2(r))dT, 

^ j=i j=i 


-Va =v;i, 
Xi2 = Ka* 


(5-11) 


First, we are mtcresled to establish a tinie-asympCotic convergence property of all particles 
towards a same position and a same velocity. For this it is convenient to introduce the centre of 
mass system, in w'hi(^ equations ( ]5.1ip become a nonautonomous system of ordinary clilfercntial 
equation, deconpled with respect to the f-th particle and with respect to the two components of 
each position and velocity vector. Then the equation of the centre of mass can be studied apart. 

We dehiie position and velocity of the centre of of a system of particles with same 
mass, as 


1sl 

Vc:m - ^ E V: 


and we consider the new variables 


X, := X. - X( M, 

V, := V. - Vcw. 

In variables ( j5.14p 15 1 the equilibrium condition becomes 

(X..V.) = (0,0), Vf = l.JV, 

moreover the following identities hold: 

jV 


Ex, = o, 
1=1 
jV 

Ev, = o. 


1=1 


(5.12) 

(5.13) 

(5-14) 

(5.15) 

(5.16) 


(5.17) 


(5.18) 
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K V(;m = (VcMi. VcMj) and Xc m = (-Vcmi,-Vcm 2 ). from qr».lip we hav-c 

v™. = Va = - I'.iW) 

I»1 tsljal 


t=i 

1 - S’- 

C(t - r)^ (NA*,i(0 - A^-Vc;M 2 (r)) dr 

= ^^(A'VcMi(i)-iVVn,a(t)) 

- ^ _('<?('-’■) (A'^.YcmiW - A'"X,:Mi(r)) dr 

= - f‘c{t-r)(NXc,ui{t]-NXcmM) dr. 

Jo 

where we have used deiimtioiib V 12^ (O^. The same holds for Vcm 2 > 

V^jiU = - f — r) - fVXc«M2(^)) 

Jo 

Ill the variables (X<s V,). taking into account ) ).'). 19^ ( S.2Q| |> equations become 

^;i = -t'c.'u. + 4 D':j>w-f’=w» 


S 


f C{i - T) V (A-,i(J) - .Y„(r) + A'( 


Ali(0 - -V(:mi(t)) dr 




S 








fi ^ 

- / C{t- r){NXa[t) - 'VXsi[7) ^ NXvu\{t) - NX(Mi{r) 

- N.Vc-ui(0 + fV'-^(:M2('r))dT 




(^'c((-T)dr^A-,i{(), 


where we ha%'e used 'i.i? j ]5.18p ^ and 

= --^CMi + + VcMi = Ki. 


(5.19) 

(5.20) 
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Similarly for aud A', 2 * 

Finally, we write the following system; 


Vri = -0Viiit)-g{t)X,i{t), 
Vr2 = -0V,2{t) - 
ixi=Vii- 

= K2, 


(5.21) 


where 


g{t)-.^N fc(t-T)dT. 

Jo 


Now we will pit)vc the uiiLTurm asymptotic stability of equilibrium ( 5.1G| | providing a suitable 
Lyapunov function for system i 5.21[ |. For simplicity, system | |5.21[ | can be written, for each 
particle aud for each component, as a planar s^'stem in the variable y = (K A): 


V = -SV-g{i)X. 
X = V. 


(5.22) 


with 


y{t):=N f C(t-r)dT. 
^0 


(5.23) 


In relation to ( ^.22) we piwe the following two theorems. 


Theorem 3. Fixed a f > 0, the system \ 5.22[ |. admtCs a Lyapunov function U{t,y) with the 
properties: 




k2\\yf<U{ty)<ki ||y||^ 


(5.24) 


b) 


for all t > i, where ki, k 2 , and ks are posi^iite constants. 
Proof. Let f > 0- we define the Lyapunov function 


(5.25) 


V{t, V,X) := {V'^ + kXV + y{t)X^)i.itl 
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whpre 


g(t) 


(5,26) 


h)ig{t) = N f C{t-T)dT, 5 := supp(f) = N / C{t-T)dr, (5,27) 

Jn t>t Jo 


0 := inf ^ = e fi, ^:=sup 0 = e 

“ d>t 

_8<£> 

— * ^ ^ 9 

i) := sup |i^| = sup-, 

t>i t>i 9 

[ 9i} 20gi}^ 

k := jnin ==,-——^-—— 

y ^ 2p^+(i/?+/3t^)2 


(5,28) 


(5,29) 


(5,30) 


Ill the following equations we consider the inequalities; 

2 “ “ 2 

Then^ because of (?(f) is an increasing function, ^(t) is nonincreasing, so 


0<O, 


(5.31) 


and finally 


_8i£>o _£il2 

= -e i -p + e 2 p 


= (i-£^ 

<0, vt > t. 


(5,32) 


To prove the second mequality in a) we consider 


V{1, V, X) < + s-Y") 


where 


A.'i:=max ^1 + 1/), ^5 + . 
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To prox'c the first part of a) we consider 

Cf(c. V.X) > + QtX'^- kjp 


2 + 




.0 


5^ 


>A.*2 




^ ££' 
2 ’ 2 


where 

k'i := min 

To prove i)) ?.'e conpider the foUowmg inequalities; 

V[U V, X) = t{V^ ^kXV + gX^) + iJ [2V^{-dV' - gX) + kV^ 
^kX{-3V - gX) + gX^ + 2i/.VV'] 

= (ti - - 3k^') XV + (^5 + 05 - A:50) 

<rrT7^ TO (-2d0 + A:0} + (fc0 - XV - kgV'X^ 

< (-2d0 + kil>) + (fc|0| + ^krk) \X\\V\ - ^P0.V2 

< (-2j?0 + A:0) + + |;>f||V|-frff0X‘ 


= {-201 +hi>) 


W- 


|V| (0 + ^ 0 ) 


2 

. 2 


50 


2 


(^ + SV : 


250 


< -2dV’ + fc0 


0 + ^ 0 ) 

2p0 


V'= - 

2 




< -^*3 lly| 


where 


hi r= min 




This (completes the proof, 

Starting from Theoi'em[^we can state 
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Theorem 4. The point (V.X) = (0,0) of the Uneari::ed system | ]r).22^ u glob{dly 

uniformly asymptotically stable with exponential rate of convergence. 


Proof. Inequalities 1)5.25^ and jfS,24i imply that U satisfies the dilfereutlal inequality 


U<~v, 

«1 


By the Gronwall's inequality, 

Then, asiug agaiu I fS.24[ f, we have 

l|yWII<(^^) 


1/2 


1/2 


-I ^ ) 


(5.33) 


Rom ( 5.33[ [ conditions for the uniform asymptotic stability arc satistied svith exponeutlal 
convergence. The proof is completed. 


Retmuiug to system i ]5.21[ i, Theorem cau be applied for each particle and for each 
component of the positiou and velocity vectors. Recalling tiansfoimatlous 1 )5.14^ ( we can 
say that all agents converge. tlme-asymptoticalJy, to the position and velocity of their centre of 
mass. 


Now, let us investigate equation i|’i.i!)|i for the motion of the centre of mass. For 

(5.34) 


c'onvetkience, we write it again for the first compoueut, the second is similar: 

VcMi = -N j C{t-r) (Xt u,i(t) - ^rMi(T)) dr. 


[ 


Xchii = VcMl- 


wlth C(t - r) defined in ( rj.lQ that, from j]5.2^. gives 


C(t-T):= 




2((-t)D 8(t-T)=0= 

Then fi'om system (we get the scalar equation 

t'cMi = 

and siinilarlv 


(5.35) 


(5.36) 


VH12 = -N 


i' 


V'cM2(^) t/s dl*. 


(5.37) 
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III order to establish the asymptotic behaviour of we prove the following 

Theorem 5. Eifudlion l |o.3G| l, with inilioi datum VcA/i(£o). has a uniqun soluhon m [£o,+' 2 o). 

Proof. The existeuec of the eoliitlou is ensured by the existence of the solution proved 
for system 1 )5.i[ i. 

Since the linearity of the equation it is enough to piove that the unique solution svith zero 
Initial datum is the constant zero solution, In fact, if we consider the case of nonzero initial 
datum and we suppose that there exist two solutions, ^(£) ft»H in correspondence of the 
same initial condition ^(£q) = 2 (£q)i their difference ~ y{t) - 40 is also a solution with 
zero initial datum. If ai(£) must be necessaiily the zero solution, then j/(f) and 2 (f) have to 
coincide, ensuring the imiqueuess. 

First, v.'e observe that C(f — r) in ( 5.35[ i Is bounded, we say 

0 < iVC(i - r) < i', 

with u positive constant. Then, from 1we write the Inequalities 

~\ dt \ 

I It 


di 


= 1 ' ^ |Vcui(a)| ds, 

Jo 


where we have changed the order of integration in the integrals. Integrating both sides from fu 
to t, changing again the order of integration, and lixed a 7" > fo- we state 


IVniiCOI 5 l'^c*Mi(fo)||V'rui(a)| d.s 

Jtu Jo 

< |VcMl(fo)| [ [ ^ |H ui(fl)| 

Jo Jo 

= \Vc\u{to)\ + I'J 5(f-^) |VcMl(3)| 

< IVcMi(fo) \ + i'T J ^1 V( m: {s)\ da, Vto < t 


dr 
ds dr 
ds 


<r. 


By Grouwall’s Letmna 


yTl' 


|Vtwi(t)|<|V’c,'Mi(fo)|fi“5". vto<«<r. 

so that, due to thearbitraiinessofT, if VcMi(£o) = 0. VcMi(0 will be the constant zero solution. 
This completes the proof. 


From Theorem O descends 

Theorem 6. Each solution of equation i |5.3fi| i tenda to zero as t —* +oo. 

Proof. If the initial datum Is V^mi( 0) = 0 the unique solution is the constant zero 
solution. We consider now Vomi{ 0) > 0. In this case the solution is greater or equal to zero in 
all (0, +eo). In fact, if there is a T > 0 such that V‘bMi(£) = 0, Theorem apphed at £o = would 
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give VfcMi(t) s 0 for all f > Thou the right hand side of i !j.3Gf( implies Vcvn < 0. so V^Mi(f) 


Is a (lecrea.sliig function on (O.+co), and there exits its limit 0 < / < +co, as f ^ +oc. Hence, 


since from ( 5.3G[ 1 l^cMi (f) ^ inonotoue and Its limit exists, we have also hmt_»+oo V^CMi(0 = 0« 
NoWs if we suppose / > 0, Vc:Mi(f) > ^ for all t > 0. and from holds 


C{t-T){t-T)dT 

Ja 

= -Nl [ zC{£)ds. Vf € [O.+co), 

^0 

where we have perforinetl the change of %‘ariables i = £ — r. Taking the limit for £ 
both sides of we would have the contradiction 


(5.38) 

(5.39) 
+00 on 


0= lim VcMiCf) <-h < 0- 

(—*+oo 


where 



So the only chance is i = 0, 

The case Vomi( 0) < 0 can be treated in a similar way. The proof is now completed. 


Theore m i s t'alid also for the secoud velocity component whose equation is 

analogous to ijO^, 

Theorems ^ aud [5 eicsuie a couditiou of time- 


Finally. returning to system Il5.ll| 


asymptotic Hocking such as stated in befiiiition[T] Moreover, we'^ve also the stronger couditiou 
that all particles converge asymptotically to their centia of mass and the velocity of the centre 
of mass decays to zero. 


C. Numerical simulations 

In this section we present some numerical simulations to show the dyuamical behaviour 
of the model introduced in Section]^ Numerical results are compared with the analytical ones 
presented in Section and some similar cases are considered. 

Nimieiical tests are performetl and shown in noudimeuslonal form using the following 
dimeusioidcss quantities; 


X’:=|, /■:= 

• 'y/max _ D 


i 


Uax 
c'_ 


i 


v' 
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2d 


where /max IS the maximum conceutratiou of signal /. With these defiiiitloiis, system | ]2.ip 
be written as 




N 






X, = V.. 


( 6 . 1 ) 


iv 


a,I = DAf + iJ2\Bix,M-f. 

where we liave dropped, for simplicity, the asterisks for the iioudimetLsioiial qiiautities. Notice 
dial, due to the choice of f? as diaracteristlc length, the dimeusiotdess particle ray turns out to 
be a unit value. 

G.l. Numerical methods. The numerical approximation sdieme used here employs a 2D finite 
ilifiercuce method on a spatial domain SI := [a, 6] x [c,d] with periodic boundary conditions. 

For the parabolic equation in order to eliminate the stiff term —/we perform the 

classic exponential transformation, and we apply a central difference scheme in space and the 
parabolic Crank-NicolHon scheme in time. 

For equations ( ]G.l^ i we adopt a one step IMEX method, putting in implicit the term 
depending on the velocities aiid in explicit the gia<lk'ut term j Hundsdorfer and Verwer[ 2003^ . 
For further details we refer to Appendix [X| 

G, 2. Numerical test s. In all t he following munerlc al tests we set N par tides in a spat ial domain 
n = (0,50j X (0,i)0j with periodic boundary conditions, and we choose a suitable time interval 
of obseri’ation [0, Tj. For the initial data we fix 


and, for ! = 1,,.., M 


/(x.0)=0. 
X,(0) = X,o. 

v,d)j =V,o. 


In particular X,i) Is chosen as a rau<lom vector, such that all the particles at f = 0 are contained 
in a suitable uutial region, hxed in the domain. Then Vju = (V^ocos^i,, V^oslu^*) chosen 
with random numbers in [0, Vo^,n^], and 6^ random uumbois in (0.2ir]. Let us now describe 
some meaningful tests. 

Test L In this test we set the parameters c = 0,5. d = 5, = 2 x 10^, £> = 2 x 10^, 4 = 0,5, 

Vi, n,.v = 3, and we cuuslder A* = 10 par tides located in aa in Figure 


fi.l 


'a). The time 

interval of observation is [0,500j. Spatial and temporal discretizations are given respectively by 
Ax = Ay = 0.25 and At = 10"'^, 

Figure [iTT] shows four time stops of the numerical simulation. Here and in the next tests 
for each time step wc plot, on the left the chemoattractant concentration /(x,f). while on the 
right the positions and the velocities of the particles in the spatial domain. The lifd square 
al t = 0 is the region In which the initial positions are taken. The red marker hidicates the 
centre of of the system, and the blue anuws are the velocity vectors. We observe an initial 
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Stage iu which the particles teud to move somewhat aligned about until f = 5 (Figure pTII (b)), 
then they hegiu to cou%‘crge to their ceutre of mass about at t = 30 (Figure (c)), tiiially all 
particles stop iu a same positiou (Figure |G.1| (d)). 

I a) is shown the spatial fluctuation around the ceutre of mass system 


Iu Figure G,2 


jV 


Fi.vd) i= 53 ||X.(f)-X(.u 


1=1 


as a function of the time. Iu Figure [ 6 ^ (b) is shown the velocity fluctuation aiuiiud the ceutre 
of mass 

S' 


Fiv(t) ;=^llv,(i)-Vc,-M 


For t > 61 Flx{t) and F’^ (t} are less than 10“Notice that the square root of Fl_\{t) and 
pi opart ioual to the standard deviations of XT(f) and V,(0 with I'eepect to the position 
and velocity of the centre of mass. Figure 1 ) 2 (c) displays the norm of the velocity of the 
centre of mass HVcmCOII '"crsus time. For I > f)2 this velocity is less than 7.8 x 10“^. The 
oscillating pattern, shown here and in the next ta‘'ts. can be attributed to the uumerieal error, 
as it is couflrmed using finer meshes for the discretization. Flgme [ 6 ^ (d) is obtained, for the 
same simulation^ asing a spatial and temporal steps lespeclively of Ai = Ay = 0.0125 and 
Af = 10“'^. For t > 8 S llVcM(f}ll ^ 2 x 10“^. With smaller values of Ax, Ay, and At, 

lower values of ||V(>M(t)|| are achlet’ed (data not shown). 


In the next two tests, starting from Test 1, we change only the two parameters related to the 
alignment and to the cbemotaxls. in order to hivestigate the competition between these two 
effects. 


Tent S. In this test we set 5 = *> = 10. while the other parameters are set as in Test 1. The 
time of observation is [0,3500]. 

Particles and chemotactic signal are not shown. In Figure |fi.3| (a) (b) we plot the spatial 
and the velocity fluctuatlous around the ceutre of mass sy.stem /7.v(f) and F7v(t). For t > 1680 


FIxit) and F7v (t) are less than 10 Figure |C.3| (c) displays the quantity ||VcM(i)ll 
function of time. For t > 2134 this value is less than 1.82 x 10“". 


as a 


Teal S. In this test we set ^ = 16. 7 = 10 and the other parameters ate set as iu Test 1, The 
time of observation is [0,3500], 

Tlie hgnres related to the particles and to the chemotactic signal are not shown, while in 
Figure [?r^ (a) (b) we plot the spatial aud the velocity fluctuatlous around the ceutre of 
system Flx{t) and Fl.\-{t), For t > 3109 we have values less than 10~^°. Figure^^(c) displays 
the quantity ||Vcm(£)II versus time. For t > 2579 this velocity is less than 1.87 x 10 “^. 


In Table we sunuuArlze the results of the first thi'ee testa. We can observe that Increasing 
values of the ratio ^ imply a decreasing rate of eoirvergence of tlie particles. 





EL DI COSTANZO AND R NATALINI 




Table 1, Spatial and velocil> fluctuations in the centre of uiass system, and 
velocity of the centre of mass for iliiferout %‘alues of the parameters Q aud f used 
in Tests 1 3. 


Test 

0 

7 

Time after which 
mnx{Flx^Fl\ } < lQ-^° 

Time after which 
||VrM(f)||< 7.8x10-2 

1 

& 

2 X 10^ 

61 

52 

2 

li- 

10 

1689 

1021 


IS 

10 

5155 

itU 


Test 4- fJi this simulation we fix the parameters as in Test 1. aud we double the number of the 


iuteractlug agents, considering = 20 cells bs iii Figure* 5,5 (a). The time interval of obseiA'atiou 
is [0,500]. 

In Figuie [lj7j| wL> can observe four diti'erent time steps of the mimexical simulation allowing 
the aggregation of the initial gioup and the convergeuce to zero of its velocity. 

In Figures [fTol (a) (b) we plot, tlie quantities F/xiO J^ud /’ir(£). Hei'e, for t > ,34 we 
have Milues less than lO"*®, 111 Figure |(b^(c) we show the quantity ||Vf;^(f)|| versus time, with 
values smaller than 8.30 x 10“’ for £ > 28. Comparing the results of this test with those of Test 
1 we can state Chat, svith the same parameters, an increasing number of cells euliances the rate 
of convergence due to the greater aiiiouiit of the expressed dieuioattractaut. 


Teal 5. In this test we consider a case iu which, for the pure Cucker-Sinale uiodel jl.l^ . i.e. 
■> = 0 iu our model, the fiocklug l)ehaviour does not occur. We fix the other parameters <j = 0,8. 
^ = 5, D = 2 X 10^, $ = 0.,5, Vb.max = 3, and X<) as in Figure [fiT] (a). The time interval of 
olkservatiun is [0,15]. Spatial aud temporal discretizations are as in i'est 1. Clearly in t.hig ease 
equations ijj^i^ aud ij^:j are decoupled. Since tj > 1/2, aud taking iuto account the initial data, 
according to the r(*siiltK iu Ha and Liu ( ]2QQ9 1 , the Hocking of the system is not guaranteed. 

Figure [fTT] hIiowh out numerical simulation at three time steps. We can olxserve, in fact, a 
dispersion of the initial group of particles. 

In the next test we will show that, adding the chemotactlc effect, we can recover the 
time-asymptotic couvergouee of the luJgratiug group. 

Tent 6. In this test we set the same paiameters and initial data of Test 5. moreover we introduce 
a non zero value for the chemotactlc effect: ■> = 10^, The inter%'al of observation, spatial and 
temporal discretizations are fixed as In Test 1. 

Figure [6^ shows the uumerical slmulatious at four different time steps. We cau uote that 
the initial group of particles does not dispei'se. but converge in position and velocity. 


In Figures 1^ (a) (b) we plot, as iu Test 1, respectively Flx{t) and F’^'(£). For £ > 107 
these quantities are less than 10“*^. Figure 6,9 (c) shows the quantity ||Vc'm( 0II versus time. 
For £ > 65 we have values smaller than 2.1 x 10“^. 


Test 7. luspired by the ruathematleal model proposed iu Di Costanzo et al 
zebrafish lateral line development, we consider a siniulatiou in which two kin 


for the 
oT cells are 


involved: the leaders that produce the diemotactlc signal and the followers that do not produce 
auy signal, both subjected to the abgnment effect and to the attraction of the diemlcal gradient. 
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For the uumptical simuIatioQ set j = 0.5s 5 = 5, = 1.5 x 10^, D = 2 x 10^, ^ = 3, 

lo,n*AJc = 0,3. iiud Xi) as iu Figure pTIOl (a). We consider a time interval of obsert'atiou [0,.'500]. 
Spatial and tenipural discretizations aie as in Test 1. 

Figuie [CTini shows four time steps of our simulation, A single leader cell is marked In green 
colour, wliile the other follower cells are in red coloin. Here the centre of mas.s is marked in 
blue. We observe, about at f = 15. that cells begin to be attracted toward the chemoattractant 
source (Figure [6.IQ (b)). In the next time steps the cells tend to converge in the centre of 
and then they stop (Figure 16.101 (c) (d)). 


In Figures G.ll (a) (b) are shosvn the spatial and velocity Hiirtuations, Fi.v(0 ^d Ffv*(t). 


aiuiind the centre of mass system. For t > 131 these quantities ate less than 10 Figure 6.11 
(c) displays the quantity ||Vc;u(0ll' ^ ^ ^28 we have x’alues smaller than 5,3 x 10 


-3 


Test 8. In this test we C'onsider a case in which the particles have zero initial velocity. If we 
C’ousidered the pure Cucker-Smalc model ) ]l.ip , we would have the particles constant in the initial 
pasltlou with zero velocity. Here we show that the chemotactic effect ensure the convergence of 
the group of particles. As in Itest 7 we consider a single leader that produce the chemoattractant, 
aud other followers that do not produce any sigual, but they foUow the gradient produced by 
the leader signal. 

For the numerical simulation we set <t = 0.0, = 2. f = 10^, D = 2 x 10^, 4 = 3, 

lo.niAJc = 0. j^ikI Xo as in Figure 6.12| (a). The time inteiwal of observation is the same as Test 
5. Spatial and temporal <iisi'rotizatians ai'e fixed ss in Test 1. 


Figme 6.12 shows four time steps of our simulation. Green colour marks the leader ceD, 
while red colour is for the followers. The centre of mass is maiked in blue. In early time steps 
wc obseiA'e the velocity vectors that point the source of chemoattractant (Figure [6.12 (b)). In 
the next time steps all particles go to com'crgc in position au<l velocity. 


In Figmes 6.13 (a) (b) we plot respectively Fl.xit) and versus time. For t > 67 we 

obtain values less thau 10“*”, In Figme 6.13 (c) is shown the quantity || Vcu(t)ll ^ ^ function 
of time. For t > 109 we have values smaller than 4.3 x 10“^. 

Tent 9. In this test we simulate the system only under the diemotactic effect, ueglectlng the 
alignment term, that is we set 5 = 0. The other paiainetcrs are *f = 10^, D = 2 x 10^. 4 = 1.5, 
and the initial data Va.mAx = 0.8, Xj as in Figure |li.i4 


[0,4000]. Spatial and temporal discretizations are fixe 


(a). The time interval of obser^mtiou is 
as iu Test 1. 


In Figure [fm| we plot fom' time steps. Although the Initial group aggregates, we do not 
oljsei VC. in our tiiuo of observation, a conveigeuce of the partirles. rather they show an osclllathig 
behaviom aioimd their centre of mass. 

In Figures [6.1.S| (a) (b) we plot respectively aud F’ii (I). In this ease, the spatial 

aud velocity duct nations around tlie centre of mass remain boimdod but do not converge to zero. 
In particular, in the time interval of observation, we ha%'c Flx{t) > 0,88 and, if we consider a 
trend Ihie, its slope seems to deciesses monutouically approaclilng to zeio. For example, in 


relation to Figme G.IS (a), the slope of the linear' til ou the values of computed ou the 

time inter^nis |0,SO0] and [3200.4000], changes from -1,46 x 10“* to —1.84 x 10“'*, Morem'er, 
we remark that in previous simulations, eoutahilng the alignment effect, ou time intervals much 
smaller than this test we have obtained %‘alues of Fl.x{t) less than 10“*^. 

In Figure (c) is sliowu |1 V(’m( 0II ^ function uf time. For t > 254 we have values 

smaller than 3.79 x 10“*. Performing the same numerical test with the finer discretization 
Aj: = Ay = 1.25 and At = 10"^, we find that, fort > 120. ]]VcM(t)|| is Bmaller than 1.41 x 10"^, 
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sec Figure (1.1.') (d). From this we cau deduce that the velocity of the eeutre of mass goes to 
zero. 

Om uumoiical results suggest that our model, in absence of alignment and with the only 
chemotaetic effect, is unable to reproduce stationary patterns. 
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(A) 





(d) 


Figure 6.1. Test 1. Numerical simulation in a spatial domain li = [0,50] x 
[0,50] with periodic boundary conditions, and in the time interval [0.500]. The 
paiainetei'K values are tr = 0.5. ^ = 5, *y = 2 x 10-. = 2 x 10^. ^ = 0.5, 

Ib.tnAx = 3. and Xo randomly taken in the red sqnaie shown in (a). Spatial and 
temporal disrietizations are respectively Aj: = Ay = 0.25 and At = 10“'*. The 
four plots arc respectively at time step.s t = 0.5.50,500. On the left there is the 
chemoattractant concentration /(x,f)t while on the right the positions and the 
%‘elocities of the particles. The red market indicates the centre of mass of the 
system, and the blue airuws are the velocity vectoiis. 
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(0 (d) 

Figure G.2. Test 1, (a) Spatial fluctuation around tlie centre of ma.s.s system 
Flx{t) := ||X,{f)-X(;M(t)||^ as a fimctiou of tlie time, (b) Velocity 

fiuctuatioTi ai'uuud tUo centre of mass systeoj FU'{t} := i II(t) — Vcm( 1) 11^ 
as a function of the time (x-axis shows only a part of the time domain). For 
t > 61 Flxit) and Fly{t) are less than lO"*®, (c) Norm of the velocity of the 
centre of mass ||Vcm( 0II ^'^rsus time (only a part of the time interval is shown 
on the x-axis}. For £ > 52 tins velocity is less than 7.8 x 10“^. The plot Is 
related to a numerical simulation with spatial and temporal steps respecti\'ely of 
Aj: = Ay = 0.25 and At = 10“'*. (d) With a hner mesh, Ax = Ay = 0.125 and 
At = 10“^. smaller values of || Vciij(0ll obtained: for £ > 88 the maximum 

value in the plot is less than 2 x 10“^. This confirm that the oscillating pattern 
shown in (c) can be attrilmted to a numerical error. 
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(b) 


(b) 



(c) 

Figure G.3. Test 2. Spatial aud velocity fluctuatiuus around the ceutre 

of ma.ss system Flx(t) aud F^'(0. For f > 1G89 are less thau 

10“^® (ou the X'lixis only a part of the time doiualu is shown), (c) Norm of the 
velocity of the centre of mass ||Vc*M{t)|| wiisUH time. For t > 2134 this velocity 
is less than 1.82 x 10"^, 
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(a) 


(b) 



(c) 

Figure G.4. Test 3. (a)-{b) Spatial aud velocity fluctuatious around the ceutre 
of mass system: Flx{t) and For t > 3109 FTv(0 and Flv{t) are less 

thau 10“*^ (ou the x^axis only a part of the time doiualu is showu). (c) Norm 
of tUi* velocity of the centia of mass |[Vt'iij(t)|| veisus time. For t > 2579 this 
velocity is less thau 1.87 x 10"^. 
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» I 


(d) 

Figure 6.5. Test 4. Numerical simulation iu a spatial domaiu li = [0,50] x 
[0,50] with periodic lioimdary c*onditions, and in the time interval [0,500]. The 
paiainetei'K values are tr = 0.5. ^ = 5, *y = 2 x 10-, = 2 x 10^, ^ = 0.5, 

Vh.cnAx = 3, and Xo randomly taken in the red sqiiaie shown in (a). Spatial and 
temporal discretizations are respectively Aa" = Ay = 0.25 and At = 10“ h The 
four plots aie respectively at time steps t = 0.10.30,500. The chemoattractant 
concentration /(x.t) Is on the left, while on the right there is the positions and 
the velocities of the particles. The red marker marks the centre of mass of the 
system, and the blue arrows are the velocity vector's. 



















35 


EL DI COSTANZO AND R NATALINI 



Tfna 



(b) 


(b) 



(c) 

Figure G.O. Test 4, (a) (h) Spatial aud velocity fluctuations around thecputrp 
of inass system Flxit) and F'iv*(0 (x-axis shows ouly a part of the time domaiu). 
For t > 34 we have values less than 10“*°. (c) Norm of the velocity of the ceutre 
of ||VcM{f)|| versus time (only a part of the time interval is shuwu on the 
x-axis). For t > 2S this velocity is less than 8.39 x 10“^. 
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Figure 6.7. Test 5. Nxmierical siuiuJation in a spatial doinaiu = [0, S0| y. 
[OiSO] with periodic boundary couditlous, and iu the time interval [0,15]. The 
parameters values are u = 0.8. ^ = 5, '7 = 0. D = 2x 10^. ^ = 0.5, Vb,io*x = 3, 
aud Xo randomly taken in the red square shown in (a). Spatial and temporal 
diKcretizatiuus aie respectively Ax = Ay = 0.2^ and At = 10“ ^ The red marker 
indicates the centre of mass of the system, and the blue arrows are the velocity 
vectors. On the left there is the chemoattractant concentration, while on the right 
the positions aud the velocitieR of the particles. In this case the motion equations 
aud the chemoattractant equation are decoupled aud for the particles we simulate 
the pure Cucker-Smale model j ]l.I[ i. Taken into account the parameters values 
aud the initial data the docking behaviom is not ensured by results m|Ha and Liu 


j2Q09[ f. In fact, from the three plots, taken r'espettively at time steps t = 0.S> IS, 
we can observe a dispersion of the initial group. 


























EL DI COSTANZO AND B. NATALINI 



Figure 6 . 8 . Test 6. Numerical simulation iu a spatial domaiu 11 = [0, S0| x 
[OiSO] with periodic boundary oouditious, and in the time interval [0,500]. The 
paiainetei'K values and the initial data ai'c as in Figure |G.7| tr = 0,8, = S, 

D = 10^. $ = 0.5, Vb.aiax = 3, and Xo as in Figure jol?] (a). Tlnu here a nouzeiu 
value for *r is diosen: 7 = 2 x 10^, Spatial and teinpuial dlscretizatious ate 
Ax = Ay = 0,25 and At = 10“'^. The four plots are respectively at time steps 
? = 0,10.30,500. On the left there is the chemoattractant coneentratiou /(x. f), 
while on the right the posltious and the velocities of the particles. The red 
marker Is the centre of mass of the system, and the blue arrows are the velocity 
vectors. Compaiing this figure with Figure |G.7[ we observe that, introdiichig the 
cheuiotactic effect, the time-asymptotic convergence of the particles Is recovered. 
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(C) 

Figure 0.9. Test 6, (a)-(b) Spatial fluctuation v'clocity fluctuation 

around the centre of mass nystcm as a fimctlou of tlie time (x-axis shows 
only a part of the time iuter\'al). For f > 107 F'iA'(f) F7k’(t) arc less than 
10“^® (ou the x-axis only a part of the time <iomalu is sUosvn). (c) Velocity of the 
centre of mass l|VcM(f)|| versas time. For t > G5 we have x'alues smaller than 

2.1 X 10-^. 
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(d) 

Figure G.IO, Test 7. Numerical siinulatlou in a spatial domain n = [0.50| x 
[0i50] with periodic bouudary cooditiorts, and m the time interval [0,50fl]. The 
paiaineters values are e = 0.5. ^ = 5, > = 1.5 x 10'^, D = 2 \ 10-. $ = 3i 
Vh.coAx = 0.3. Aiid X<( randomly taken in the red squaie shown in (a). Spatial 
and temporal discretizations are respectively Ax = Ay = 0.25 and At = 10“'^. 
The green cell marks a leader cell (•}, that produce the chemical signal, while the 
other red cells are the followers (•}, that do not produce any signal and follow the 
chemoattractant gradient. The blue marker is the centre of mass of the system, 
and the bine arrows are the velocity vectors. The fom' plots ai'e respectively 
at time steps t = Q. 15,40.500, On the left there is the chemoattractant 
concentration, while on the right the positions and the velocities of the particles. 
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(c) 


Figure 6.11. Test 7. (a)-(b) Spatial aud velocity fluctuations, Flxit) aud 
F7i'(f), around the centre of mass system (only a part of the time domain is 
shown ou the x-axis). For t > 131 and F7y(f) are less than 10“*^. (c) 

Velocity of tlie centre of ||Vcm( 0|| aa a function of the time. For t > 328 
we have values smaller than 5.3 x 1Q~^. 
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Figure G.12, Test 8. Numerical simulatlou iu a spatial domain H = [0.50| x 
[0i SO] with periodic boundary conditions, and iu the time iuteiA’al (0, SOOj. Spatial 
and temporal discretizations are respectively Ax = Ay = 0.25 and At = 10“'^. 
Iu this case we fix a zer'o initial velocity, Vo.mAx = 0- itud the parameters cr = O.G, 
5 = 2. > = 10'^. D = 2 X 10^, € = 3. and randomly taken in the red square 
shown in (a). As in Figure Q.IO we consider a single leader cell {•) and other red 
follower cells (•), Because of the zero initial wlocity, if we considered the pure 
Cucker-Smale we would have all particles constaut in time iu the initial 

position. Conversely, iu the fom' plots at time steps t = 0,5,30.500 we observe 
the time-asymptotic convergence of the migrating gi'oup. 
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(b) 



(c) 

Figure G.13. Test 8, (&)-(b) Spatliil aud velocity fiiictuatioiiH around ti»e cputrp 
of mas system as a function of the time: F’^v(^) «Did F/v (0* ^ ^ F/.v(0 

and Fl\ (t) are less than 10“(c) Velocity of the centre of mass ||Vc;u(0ll 
versus time. For t > 109 this quantity Is smaller than 4.3 x 10“-. 











































































































































































EL DI COSTANZO AND B. NATALINI 



(d) 


Figure G.14, Test 9. Numerical simulation m a spatial domain 12 = [0.50] x 
[0i SO] with periodic boundary couditious, aud iu the time iuter^’a! [0.4000] (plots 
ou the light shown only a part of the spatial domain). Spatial aud temporal 
dibcretizatioua aie respectively Ax = Ay = 0.25 aud At = 10~*. In this test 
only the chemotactic force is considered, neglecting the alignment effect (^ = 0). 
For the other values we fix y = 10^, D = 2 x 10^, ^ = 1.5, Vo^nuix = 0.8. and Xo 
randomly taken in the red square shown in (a). The plots taken at time steps 
/ = 0,10.2000,4000 display the aggregation of the Initial group of particles. Iu 
this ease the couvergent'e is not observed, but an oscillating lK*Uavioui ai'uimd 
their centre of mass is shown. 
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FiGi’RE 0.15. Test 9, (a)-(b} Spatial aud velocity fluctuations, Flx{t) aud 
F’/r(J)- aroiiiiri the centre of mass system versus time* (only a pait of the time 
domain is shown ou the X'Sxis). lu this case F’ix(i) Aiid Fl.y{t) remain bounded 
but do not converge to zero. In particular wc have Flx{t) > 0.8S aud. if vw» 
consider a trend line, its slope seems to decrease mouotonically. In particular in 
the time intervals [0. flOOj aud [3200,4000]. the slope of the lit line diauges from 
-1.46 X 10“^ to -1.84 X 10“'*, (c) Velocity of the centi-e of mass ||V(-M(t)|| versus 
time, using a spatial and a temporal discretization given by Aa* = Ay = 2.5 and 
At = 10“** (x-axis shows a pait of the time interval}. For t > 254 this quantity 
is smaller thau 3.70 x 10“^, (d) With a finer mesh, Ax = Ay = 0.125 and 
At = 10“'^, we obtain smaller t'alucs: for t > 120 we have \'alues less thau 
1.41 X 10-^, 
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7. Conclusions 


Id this paper we have proposed an extensiou of the Ciu^r-Smale model introducing a 
mathematical model for collective motion dri\'eu by an alignment and a chemotaxls effect. We 
have adopted a hybrid description, discrete for the particles and for the motion equations, and 
coutuinoiifl for the molecular level, containing the equation for the chemical sigtial 

We have studied om model by botli an analytical and a niunerical point of view. By 
the analytical point of view, using a fixed point theorem, we have proved local existence and 
uniqueness of the solution of the nonlinear system. Then, through a principle of continuation 
of solutions, we have extended such result proving global existeuce au<l uniqueness for all 
times. Moreover, we have investigated the asympCotk behaviour of particular equilibrliiiu 
configurations, corresponding to a state iu whidi all particles are located in a same position 
with zero velocity. On the linearized system aioimd such configurations we have proved, using a 
Lyapunov funct ional approach, t he asymptot ic convergence of the part Icdes in their cent re of 
with same velocity. Then the velocity of the centre of muss Is proved to go time-asymptoticully to 
zero. By u numerical point of \Tew this property has been tested ou the nonlinear system, fiiiding 
a complete coucordauce with the analytical results, aiid some numerical simulatlous have been 
proposed. Among tliese, we have discussed the couipetltlou between alignment aud diemotactlc 
effects varying the parumeters of the system. We have found a decay in the rate of convergence of 
the particles when the strength of the alignment term mcrea.ses with respect to the chemotaxis. 
On the other hand, a growth in the rate of convergence can be found with the same parameters 
aud with a greater number of interacting particles. Then we have considered the presence of two 
cell populations, inspired from the model for the zebrafish lateral line proposed In Dl Custanzo 
et al j ]2ni5^ . finding the convergence of cells towards the sources of chemical signal. Finally, 
the case of a pure chemotactlc effect, without alignment, has been simulated. Nrrmericul results 
have shown the absence of cou^’crgenee iu the migratiug group, suggest mg that our model, with 
the only diemotaxis mechanism, is not able to destribe biological procesiifts leading to stationary 
aggregates, 

Firtiue perspectives can concern the exterrsion of the stability result, proved on the 
linearized system, to the full nonlinear model. Furthermore, other iutcractions could be takou 
luto account irr the collecti%'e motiou, such as adhesion-repulslou, damplug or lateral Inhibltiou 
terms, similar to those introduced in the aforesaid paper ou the zebrafish primordium. This 
would be interesting in view of studyiug, in an analytical framework, the morphogenetic process 
arising iu tire lateral line development, and establish results in relatlou to the neuromast 
formation aud deposltlou. 


Appendix A. Numerical approximation 

The numerical approximation sdieme used in the dynamical tests employs a 2D finite 
difference method ou a spatial domain 0 ;= [ 0 , 6 ] x [c, li] with periodic bormdary coudltlorrs. If 
Ax and Ay are the spatial steps and At the time step, wc define the grid points (im; Vni fA‘)i 
where Xm = riiAx. Xn = riAy and = kAt. The notation will denote the approximation 
of a function u(a:,y,f) at the grid point 

For the parabolic equation j )6.1[ fj. in order to eliminate the stiff term —/we perform the 
exponential transformation 

f{xA) = s"^u{x,i). 
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where ii(x, t) satisfies the equation 


N 


dtu = DAu + ^ \B(X^.L)' 


(A.l) 


with initial condition 


u(x.O) = 0. 

and periodic boundary conditions. 

Now, for we apply a central clilfptence scheme in space and the parabolic C’rank- 

Nicolson scheme in time. Namely we write 

= J ^ {oy + ny 




where the second finite differences D^u and D^u ai'e given by 

j-iS. fc ~ + ^»n+l,n 

L/jil ,— 


r>2..fc . . - 2Um,n + 

- ^ — 


(A,2) 


i,n+l 


and 


_ / 1. if €B(X*.l); 

XBiX^A) \ 0 . otherwise. 


For equations ^0.1^ i we adopt a one step IMEX method, putting in implicit the term depending 
on the velocities and in expficil the gradient term ( ]Hundsdotfer and Verwei] 2003 [ i. The scheme 
rends: 


fk*[ 




At 

Xfc+i_x* 


(A.3) 


At 




where 


( f* - f* f* - f* \ 

Jm+l,n Jm.n Jm.n I 

2Ai ’ 2Ai; }' 


and represents an interpolation of the \‘aIuos Vm,n/^ oa the grid points nearest 

to Xf. ' 

About the spatial and temporal steps employed in the simulations we chose A.r = Ay and, 


in the fiiut time steps, the paraboUc CFL condition At 


Ax^ 

IT 


, which has sbowu more stability 
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lu the numerical tests, taken int-o account in | )A.ip the discontinuous source term aud the value 
of u{x, f) initially near to zero. In practice, after few time steps, we have release the parabolic 
CFL in favour of At ~ Ax. 
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